function T = rctrap(f,a,b,n)
% T = rctrap(f,a,b,n)
% Quadrature using the recursive trapezoidal rule.
% f is the name of the function, input.
% a is the left  endpoint, input.
% b is the right endpoint, input.
% n is the number of times for recursion, input.
% T is the recursive trapezoidal rule list, output.
m = 1;
h = b - a;
T = zeros(1,n+1);
T(1) = h*(feval(f,a) + feval(f,b))/2;
for j = 1:n,
  m = 2*m;
  h = h/2;
  s = 0;
  for k=1:m/2,
    x = a + h*(2*k-1);
    s = s + feval(f,x);
  end
  T(j+1) = T(j)/2 + h*s;
end
